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■ Abstract: The rate for the fusion process p + p d + e+ + i^e is calculated using non- 

, relativistic effective field theory. Including the four-nucleon derivative interaction, results are 

; I ' obtained in next-to-leading order in the momentum expansion. This reproduces the effects of the 

effective range parameter. Coulomb interactions between the incoming protons are included non- 
perturbatively in a systematic way. The resulting fusion rate is independent of specific models 
' and wavefunctions for the interacting nucleons. At this order in the effective Lagrangian there 

is an unknown counterterm which limits the numerical accuracy of the calculated rate given by 
the squared reduced matrix element A^(0) = 7.37. Assuming the counterterm to have a natural 
^ \ magnitude, we estimate the accuracy of this result to be 6% - 8%. This is consistent with previous 

nuclear physics calculations based on effective range theory and inclusion of axial two-body weak 
. currents. The true magnitude of the counterterm can be determined from a precise measurement 

\ of the cross-section for low-energy neutrino scattering on deuterons. 

o 
o 

1 Introduction 

V* ' 

^ ' One of the most important problems in modern physics is the nature and properties of 

C ■ neutrinos. These were for a long time thought to be massless and stable, but experiments 

during the last decade have consistently shown this to be incompatible with the observed 
neutrino oscillations [ffl]. Historically and even today the fusion processes in the Sun are 
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^ • among the few available and abundant sources of low-energy neutrinos available for exper- 



imental investigations. In order to study oscillations in the detected fluxes, one needs to 
be sure of the production rates in the different nuclear reactions taking place in the Sun. 

The basic process is proton-proton fusion p + p ^ d + + I'e- It was explained more 
than sixty years ago by Bethe and Critchfield when nuclear physics was still at a very prim- 
itive stage[^]. When the field had matured, it was reconsidered in the light of more modern 
developments by Salpeter who included effective range corrections [^. Applications to the 
specific conditions we have in the Sun were investigated by Bahcall and May|0]. This work 
was later extended by Kamionkowski and Bahcall who also included the effects of vacuum 
polarization in the Coulomb interaction between the incoming protons||5|. In spite of the 
enormous progress in nuclear physics during this time, the methods and approximations 
made in these different calculations were essentially the same with a resulting accuracy 
in the fusion rate of a few percent. Including strong corrections due to mesonic currents 
at smaller scales, the uncertainty in the rate is now around one percent H]. This is very 
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impressive for a strongly interacting process at low energies where ordinary perturbation 
theory cannot be used. 



In the light of the importance this fundamental process plays in connection with the 
solar neutrino production and possible neutrino oscillations, it is natural to reconsider 
the process from the point of view of modern quantum field theory instead of the old 
potential models used previously. A first attempt in this direction was made by Ivanov et 
al.[^. In their relativistic model they obtained a result which was significantly different 
from the standard result based upon conventional nuclear physics models. Subsequently 
it was pointed out by Bahcall and Kamionowski that their effective nuclear interaction 
was not consistent with what is known about proton-proton scattering at low energies 
where Coulomb effects are important |8|. In a more recent contribution this defect of their 
calculation was removed and better agreement with standard results have been obtained|^]. 

The approach of Ivanov et al. is based upon relativistic field theory and should in 
principle yield reliable results. But it is well known that it is very difficult to use con- 
sistently a relativistic formulation for bound states like the deuteron. In addition, the 
fusion process considered here takes place at low energies and should therefore instead be 
described within a non-relativistic framework. Then all the large-momentum degrees of 
freedom are integrated out and one is left with an effective theory involving only the phys- 
ically important field variables. The underlying, relativistic interactions are replaced by 
non-renormalizable local interactions with coupling constants which must be determined 
from experiments at low energies. Along these lines the proton-proton fusion rate has 
been calculated by Park, Kubodera, Min and Rho using chiral perturbation theory in the 
low-energy limit []lO| . They obtain results in very good agreement with previous nuclear 
physics calculations. This is to be expected since they make use of phenomenological nu- 
cleon wavefunctions which fit low-energy scattering data very well. The drawback is that 
the results cannot be derived in an entirely analytical way. 

A more fundamental approach to nucleon-nucleon interactions at low energies has been 
formulated by Kaplan, Savage and Wise in terms of an effective theory for non-relativistic 
nucleons[[ll| [|12|. It involves a few basic coupling constants which have been determined 
from nucleon scattering data at low energies. With no more free parameters to fit it can 
then be used to make predictions for a large number of other experimentally accessible 



quantities ||13|]. The effects of pions can be included using the established counting rules and 
higher order corrections can be derived in a systematic way. When the energy is sufficiently 
low as for the fusion process considered here, the effects of pions can be integrated out and 
absorbed into the coupling constants of the contact interactions. The resulting effective 
field theory which is sometimes called EFT('^) then involves only nucleon fields 1 14]. In 



proton-proton scattering at low energies the Coulomb repulsion has a dominant role and 



can naturally be incorporated into this theory] 15 1. As a direct result one can derive the 
difference between the strong scattering length which should be approximately the same as 
in proton-neutron scattering, and the observed one which is modified by Coulomb effects. 



The relation is very similar to the old result by Jackson and Blatt|16]. Corrections due 
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to effective-range interactions can also be included but with more difficulty due to highly 
divergent integrals involving Coulomb wavefunctions |l7|] . 

With this understanding of low-energy proton-proton elastic scattering, one can cal- 
culate the leading order result for the fusion process p + p^d + e'^ + i'e taking place at 
essential zero initial kinetic energy [^]. With the use of a non-standard representation of 
the Coulomb propagator, the result is in full agreement with the corresponding leading 
order nuclear physics result and depends only on the physical proton-proton scattering 
length. To next order in the effective field theory expansion, one can derive higher or- 
der corrections 1 19] which also have the same structure as the corresponding effective-range 



corrections from more standard nuclear physics@[^. However, at this order there appears 
an unknown counterterm in the effective Lagrangian which will enter as a correspondingly 
unknown term in the result for the fusion rate. It can be determined from other related 
reactions. The most promising is neutrino scattering on deuterons at low energies which 



has been investigated in the same effective theory by Butler and Chen|2C|. When this is 



done, we will also have a more accurate and predictive result for the fusion rate. 

In the next section we present the theoretical framework which in the following will 
be used to calculate the proton-proton fusion rate in next-to-leading order in the mo- 
mentum expansion of the effective theory. This is done both for the proton-neutron and 
proton-proton sectors. A short summary of the leading order calculation is given in section 
3 followed by a more detailed calculation of the next-to-leading order corrections. The 
derivation of the rate is completed by the inclusion of the effects of the counterterm. In the 
last section the obtained result is discussed and compared with what is obtained by other 
methods. Using a recently improved matching of the coupling constants in the proton- 
neutron sector|2^], we obtain a final result of A^(0) = 7.37 for the squared reduced matrix 
element giving the fusion rate. This is to be compared with the result of A^(0) = 7.08 
derived within the corresponding effective range approximation of nuclear physics. From 
the dependence of the result on the magnitude of the unknown counterterm, we estimate 
the uncertainty to be 6% - 8%. This is significantly more than in other approaches where 
the unknown counterterm is replaced by definite mesonic contributions with a resulting 
accuracy claimed to be close to 1%. In the present effective theory a corresponding ac- 
curacy can only be hoped for when the counterterm is accurately determined in some 
other process. Finally, in an appendix we present a new and simpler method to regularize 
divergent integrals involving derivatives of the Coulomb wavefunctions at the origin. 



2 Theoretical framework 



In the fusion reaction p + p — > d + + z^e at low energies the incoming protons are in 
an antisymmetric spin singlet state. The deuteron d has spin S = 1 and the process is 
thus a Gamow- Teller transition mediated by the weak axial current operator A_ which 
also lowers the isospin by a unit. For a given kinetic energy E and relative velocty Vrei 
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of the protons in the initial state, the reaction cross-section is then given by the standard 
formula 



a{E) = ^^f{E)mA.\pp)\^ (1) 

where the squared matrix element must be spin averaged. Here Ga is the weak axial 
vector coupling constant, rrif, is the electron mass and f{E) is the Fermi function resulting 



from the integration over the available phase space of the leptons in the final state 1 22]. 
The available energy in the process for fusion at rest is set by the neutron-proton mass 
difference AM = Mn — Mp which is 1.294 MeV and the deuteron binding energy B = 
2.225 MeV. This gives an energy release of 0.931 MeV carried away by the leptons. The 
temperature in the core of the Sun is approximately 15 x 10^ K which corresponds to an 
average proton momentum around p = 1.5 MeV and a much smaller energy. We will 
therefore in the following assume that the initial proton energy E approaches zero. The 
kinetic energy of the lepton pair will be much smaller than the momentum 7 = V BM of 
the bound nucleons with reduced mass M/2 in the deuteron. With the above value for 
the binding energy it follows that 7 = 45.7 MeV and thus to a very good approximation 
one can just ignore the momentum transfer between the leptons and the nucleons. 

The difficult part of calculating the fusion cross-section (|l[) lies in the hadronic ma- 
trix element Tfi{p) = {d\A-\pp) which is a function of the initial proton momentum 
p = \/EM. Its magnitude can easily be estimated [22]. When the proton momentum goes 
to zero, the pp wavefunction ipp{r) becomes constant over the range of the deuteron. It is 
simply given by the the Sommerfeld factor Cjj = e~'^'^^'^\T{l + ir])\ where rj = aM/2p char- 



acterizes the strength of the Coulomb repulsion between the protons 1 23]. The probability 
to find the two protons at the same point is therefore 



2 _ 27rr? 

At very low energies when r] gets large, it becomes exponentially small and is the dominant 
effect in the fusion reaction. Similarly, in lowest order the deuteron wavefunction is simply 



2tt r 

A rough estimate for the nuclear matrix element is then 



(3) 



Tji{p) = J d\Mr)M^) = C,^^ (4) 

This result sets the scale for the fusion rate. It is therefore natural to define the reduced 
matrix element 10] 
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which contains all the interesting and important physics. It is expected to have a value of 
the order of one and is now the conventional way of presenting theoretical results for the 
fusion rate. The goal of the present paper is to calculate this number in a more model- 
independent way by purely analytical methods without using any other phenomenological 
input than the scattering lengths and effective ranges appearing in nucleon-nucleon scat- 
tering. 



2.1 KSW effective field theory 

During the last couple of years much progress has been made in understanding the low- 
energy properties of few-nucleon systems from the non-relativistic effective field theory 



proposed by Kaplan, Savage and Wise|ll]. At these length scales the proton and neu- 
tron are considered to be structureless point particles described by a nucleon isodoublet 
A^"^ = {p, n) Schrodinger field. For energies well below the pion mass mjr, all interactions 
including those due to pion exchanges, will be local. In the effective Lagrangian they can 
be thus represented by terms involving only the nucleon field and derivatives thereof in 
such a way that all the symmetries obeyed by the strong interactions are preserved. At 
the lowest energies only S'-waves will contribute. Including no more than terms of dimen- 
sion eight in the derivative expansion, there are only two possible interaction terms in the 
Lagrangian parametrized by the coupling constants Cq and C2. It can be written as 

+ \C2{{N^V^IIN)-{N^II.N)^ + h.c] (6) 

where the operator ^ = — ^)/2. The projection operators Ilj enforce the correct 
spin and isospin quantum numbers in the channels under investigation. More specif- 
ically, for spin-singlet interactions 11, = (T2T2rj/\/8 while for spin-triplet interactions 
Ilj = a20iT2l \f^. This theory is now valid below an upper momentum A which will 
be the physical cutoff when the theory is regularized that way. Since the pion field is 
integrated out, all its effects are soaked up in the two coupling constants Co and C2. 
Then the value of the cutoff A will be set by the pion mass rriT^. In this momentum 
range all the main properties of few-nucleon systems are now in principle given by the 
above Lagrangian. More accurate results will follow from higher order operators in this 
field-theoretic description |]14|. 



The effective Lagrangian @ is non-renormalizable and divergent loop integrals must 
be regularized. For this purpose one can use the OS scheme of Mehen and Stewart [p^] 



which is a generalization of the original proposal by Gegelia|25]. An equivalent method is 
the PDS scheme which was invented by Kaplan, Savage and Wise[ll]|12| and is based on 
dimensional regularization. We will use it here. The a 'priori unknown coupling constants 
Co and C2 can be determined in terms of experimental quantities measured in low-energy 
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nuclear reactions. The size of the dimension-six couphng constant Cq will then be deter- 
mined by the scattering length a in nucleon-nucleon scattering while C2 is found to be 
proportional to the effective range parameter tq. Both of these coupling constants will 
depend on the renormalization mass fi which enters in the PDS regularization scheme. It 
can be chosen freely in the interval 7 < < m-,^ but physical results obtained from the 
effctive theory should be independent of its precise value. 



2.2 Proton-neutron interactions and the deuteron 



The deuteron will appear as a bound state in proton-neutron scattering in the spin-triplet 
channel. It is then natural to determine the corresponding coupling constants by matching 
the results to properties at the deuteron pole of the scattering amplitude. The residue 
at the pole gives the renormalization constant Z of the deuteron interpolating field which 
replaces the wavefunction of the bound state[12|. In lowest order of perturbation theory 



Figure 1: Feynman diagram representing the leading order contribution to the deuteron wave- 
function renormalization. The crosses represent the coupling to the interpolating field for the 
deuteron while the lines are nucleon propagators. 

one finds the renormalization constant from the irreducible 2-point function shown in 
Fig. 1. At the two vertices the interpolating field acts with energy E and zero momentum 
and a strength which we choose to be —1. In the intermediate state there is a neutron and 
a proton which propagate with relative momentum k. Integrating over all these momenta 
we then find the value of the diagram, 



(2-Kf E -\s?/M + ie 



This divergent integral is now made finite using the PDS regularization scheme[ll| and 
gives 

So(S) = -^(/i-^/^M^) (7) 
The renormalization constant 



which is evaluated at the deuteron pole where the binding energy is = thus 
takes the value Zq = y^Svi^/M at this order of perturbation theory. It is independent of 
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the coupling constant Cq whose effect must be summed to all orders in order to find the 
non-perturbatively bound state in this channel. When it takes the special renormalized 
value 

^«(^) = St^ 

we see that the reducible chain of bubble interactions in Fig. 2 just gives the same result 
as for the single bubble in the irreducible diagram in Fig. 1. In this particular channel 

<x>-o 

Figure 2: Chain of proton-neutron interactions mediated by the leading order contact term. 

one shall therefore not sum such chains of bubble diagrams when one describes the bound 
state deuteron by an interpolating quantum field. 

While the coupling constant Cq must be treated non-perturbatively, the effects of the 
derivative coupling C2 are included only to first order. The corresponding renormalized 
coupling constant is found to be CK/^) = pdM {Cq (fi))'^ /Stt where pd = 1.76 fm is the 
spin-triplet pn effective range scattering parameter evaluated at the deuteron pole[^][|l^]. 
It will also contribute to the renormalization constant Z via the perturbative diagram in 
Fig. 3 for the 2-point function. It has the value T,2{E) = C2MET,q{E) which gives the 
total contribution 

-1/2 



Zo — Zn 



1 



7M 



(10) 



In the limit // ^ 7 the dependence on the regularization mass p is seen to go away. The 



previous value Zq then gets modified by the factor v^d = l/\/l — JPd- This corresponds 
to a change of the normalization of the deuteron wavefunction (^) which now becomes 
''Pdi'f) ~^ \/^V'rf(^) ill agreement with effective range theory in nuclear physics [ 26 1. Here 
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Figure 3: Effective-range correction to the deuteron wavefunction renormalization constant. 

it is only valid at large distances since properties of the deuteron at scales less than l/niT^ 
are not accessible in this theory. Also, it is strictly only valid to first order in an expansion 
in powers of jpd since it is obtained perturbatively in the coupling constant C2 . Since 
the expansion parameter has the rather large value jpd = 0.41, it is desirable to improve 
the convergence of perturbation theory in this coupling constant. This has recently been 



achieved by Phillips, Rupak and Savage|21| whose method we will apply at the end of the 
more conventional approach we present first. 
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2.3 Coulomb interactions and the proton-proton wavefunction 



In the absence of strong interactions, the incoming proton-proton state with center-of-mass 



momentum p is given by the Coulomb wavefunction 1 27] 



oo 

Here p = pr and ai = argr(l + £ + irf) is the Coulomb phaseshift. At low energies only 
the S'-wave will contribute. It is given in terms of the Kummer function Af(a, 6; z) as 

Fo(p) = Cr,pe-'fM{l - irj, 2; 2ip) (12) 

which is a confluent hyper geometric function. 

The strong interactions between the protons can now be included using the same KSW 
Lagrangian (P) but now with coupling constants Cg and which also get renormalized [p!7|] . 
As in the proton-neutron channel one must again consider the coupling Cg to all orders in 




Figure 4: Elastic scattering due to chain of bubble diagrams with Coulomb interactions. Incoming 
and outgoing particles are in Coulomb eigenstates. 

perturbation theory. In this way one finds that proton-proton elastic scattering is given 
by the infinite sum of all chains of Coulomb-dressed bubble diagrams as shown in Fig. 4. 
Each bubble is given by the Coulomb propagator 

Gc{E;r ,r) = M ——— 

It satisfies the Lippmann-Schwinger equation Gc = Go + Gi^VcGc where Vc is the 
Coulomb potential and 

Go(i?;q',q) = , (27r)3j(q^-q) (14) 

is the free propagator in momentum space. Iterating this functional equation we see from 
Fig. 5 that it corresponds to the exchange of zero, one, two and more static photons. 
Since a single bubble in Fig. 4 corresponds to the propagation of the proton pair with 
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energy E = /M from zero separation and back to zero separation, it has the value 
Jq{p) = Gc{E-v' = {),v = [)) or 

Uv) = m[4^, '7!^\ . \ . (15) 
The integral is seen to be ultraviolet divergent, but can be regularized in the PDS scheme in 




+ • • • 



Figure 5: The Coulomb propagator can be represented by an infinite sum of exchanged static 
photons between the two charged particles. 



d = 3 — e dimensions. When contributions from poles in d ■ 
one finds lITsi 



2 dimensions are subtracted, 



Hp) 



aM2 
A-K 



7 + ^^^ + ^-2^- 



47r 



Here Ce = 0.5772 is Euler's constant and the function 



H{ri) = tpiir]) + 



1 

2ir] 



ln{iri) 



(16) 



(17) 



is known to appear in these Coulomb scattering problems The divergent 1/e piece will 
be absorbed in counterterms representing electromagnetic interactions at shorter scales. 
This replaces the bare coupling constant Cq with the renormalized value Cg (/i). It can be 
found by matching the calculated proton-proton scattering amplitude to the experimental 
one. This is usually given by the measured scattering length Up = —7.82 fm when the 
proton momentum p — > 0. Thus one finds [ 15 1 



1 



M 



+ Jo(0) 



(18) 



Since the function ( p!7| ) is dominated by its real part /i(r/) = l/{12'rf') + 0{ri^^) which goes 
to zero when p — > 0, we see from the form of Jo(0) in (|16|) that it is natural to introduce 
the ^-dependent scattering length 



1 



a(/i) 



— +aM 



1 , 1 3 



(19) 



where a is the fine-structure constant. It corresponds to the Jackson-Blatt relation be- 
tween the strong and Coulomb-modified proton-proton scattering lengths p^. Then we 
can write 

47r 1 



M l/a{ii) 



(20) 
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which is now on the same form as (P) for the bound-state case. 

In next order of the momentum expansion the derivative couphng C2 in is intro- 
duced perturbatively. Again matching to low-energy proton-proton scattering, one finds 
CK//) = pdM(Co (/x))^/87r where pp = 2.79 fm is the the proton-proton effective range 
parameter. It is not affected by Coulomb corrections to this order in the effective theory. 
However, the Cf coupling gives an important contribution to the scattering length ( |l9| ) 
which picks up an additional term —ppp/2 in the parenthesis jl7||. 

2.4 Gamow- Teller transition operators 

The dominant weak transition matrix elements in the basic fusion rate formula (||) are 
due to the elementary isovector axial current operator 

aL^^ = N^ar^N (21) 

which converts an incoming proton into a neutron with the proper spin and isospin quan- 
tum numbers. This is the ordinary one-body interaction depicted in Fig. 6a. But when 
we include the dimension-eight derivative operator in (^ higher dimension weak transi- 
tion operators must also be considered. These were first discussed by Butler and Chen in 
connection with elastic and inelastic scattering of neutrinos on deuterons ||20t| . In our case 
there is only one such operator which can be written as 

aL^^ = LiA{N^IlN)\N'^Ii^N) (22) 

where the projection operator 11 = cj2<tt2/\/8 acts in the spin-triplet final proton-neutron 
state while n_ = a2T2T-l \/8 acts on the spin-singlet proton-proton initial state. The weak 
axial vector coupling constant has been factored out so that the effective coupling constant 




Figure 6: One-body interaction in (a) represents the weak axial current vertex while the two-body 
interaction in (b) represents the higher order axial vector counterterm. 

is just L\A. This new two-body operator represents weak transitions taking place at shorter 
length scales than considered in the effective theory and the corresponding vertex is shown 
in Fig. 6b. Typically it represents transitions due to pion interactions and other two-body 
interactions. In the effective theory it will act as a counterterm which can absorb the 
dependence on the renormalization mass /i. Its actual magnitude is presently unknown. It 
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can be estimated from dimensional arguments combined with the renormahzation group. 
Even better would be to determine it in some other weak process where its contribution 
could be isolated and measured. 



3 Hadronic matrix elements 



We are now in position to calculate the hadronic matrix elements Tfi{p) = {d\A-\pp) 
of the weak transition operator. The initial proton state is constructed in terms of the 
Coulomb wavefunctions as in the elastic scattering case. For the final state deuteron we 
use the interpolating field which contains a proton-neutron state with the amplitude Z 
which is the wavefunction renormahzation constant. We will initially consider only the 
action of the axial current operator ( ^T|) . The effects of spin has been separated out and 
will not enter the following calculation. 

3.1 Leading order result 

In lowest order of the effective theory only the dimension-six operators will contribute with 
coupling constants Cq and Cg in the deuteron and proton-proton sectors respectively. The 




Figure 7: Feynman diagrams contributing to proton-proton fusion in leading order. 

transition matrix element Tjj then gets contributions from three classes of diagrams shown 
in Fig. 7. After being hit by the weak current, the proton-proton system is transformed 
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into a bound deuteron. The value of the simplest diagram in Fig. 7a is then seen to be 
ZqAq{p) where Zq is the constant derived in the previous section and 

, , , f d^k M , , 

There is a factor (-1) from the deuteron vertex and the bound proton-neutron propagator 
is — M/(k^ + 7^). In addition, we have introduced the Fourier transform ?/'p(k) of the 
Coulomb wavefunction ( pT[ ) when the protons have the center-of-mass momentum p. In- 
cluding next the strong interaction once between the two protons as shown in Fig. 7b, we 
get the contribution ZoCo-Bo(p)V'p(0) where 

is a convergent integral and the last factor ipp{0) = Cr^e^''^ gives the amplitude for the two 
incoming protons to meet at the first vertex. Going to higher orders in the coupling Cq 
we will add in Coulomb-dressed bubble diagrams as in Fig.7c. Each bubble is of the same 
form as in proton-proton scattering in Fig. 4 where the contribution from each bubble 
is given by Jq{p) in (|l5|). Adding up these diagrams, they are seen to form a geometric 
series with the sum Cg/(1 — CqJq). The total contribution from all the three classes of 
diagrams thus gives the lowest order transition amplitude Tji (p) = Zq Tq (p) where 



nip) 



(25) 



The term involving Cq can now be expressed in terms of the proton-proton scattering 
length Op in (18) and is independent of the renormalization scale ^. 



For the explicit evaluation of this matrix element it is necessary to introduce the 
Coulomb wavefunction (12). Since the first term of the momentum integral is the product 



of two Fourier transformed functions, we find that it simplifies in coordinate space to 

Ao{p) = MCris''"' / drre-'^"'+'P^''M{l-ir],2;2ipr) (26) 





MC.e^'^o ^ / 2ip 
;2Fi{l-irj,2;2; 



(7 -I- ipY \ ^ + ip 

Now the hypergeometric function 2-^1 (a, b] z) = [1 — z)'"" so that the final result can be 
written as 



Ao (p) = C^e''^' ^j^e'" (27) 
p + 7 



In the expression (|24| ) for Bq{p) we notice that the integral over k' gives the complex 
conjugate value of the Coulomb wavefunction at the origin. It therefore takes the form 

^°(^) =^]l2^'J (2vr)3k^+7^p2-q2+.e ^.(°^ 
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The integral over k is just the previous result for Aq^q) so that 



When the momentum of the incoming proton is non-zero it yields in general a complex 
result. 



In the fusion limit p — > we now find that the first term (27) simplifies to 

Ao{p^O) = C^^e^+"'° (29) 
7 

where the parameter x = otM/^. Similarly, the second term Bq{p) becomes proportional 
to the integral 

roo farctan(^) 

I{X) = / dx (30) 
JO — 1 + TT X 

in the same limit when we use x = 2TTrj{q) as a new integration variable. Repeating this 
calculation with a different representation of the Coulomb Green's function, it can be 
shown that the integral takes the value [ 18 1 

I{x) = --e^Ei{x) (31) 
X 



when expressed in terms of the exponential integral 

/■OO 

Elix) = / dt 



With Zq = y/8iTj/M for the renormalization constant, we thus find for the full matrix 
element the result 

The reduced matrix element in leading order is therefore 

Ao(0) = - aMup /(x) (33) 

This is also the canonical result from standard nuclear physics[^. The parameter x = 0.15 
and thus the integral /(0.15) = 4.96. Combined with the measured value Op = — 7.82fm 
for the scattering length, we then have Ao(0) = 2.51 for the reduced matrix element. In 
the formula for the fusion rate it gives the contribution Aq(0) = 6.30. From previous 
applications of the effective theory |jl^, we know that leading-order results are typically 
within 20 - 30% of the correct values. Going to next order in perturbation theory, the 
accuracy is expected to increase to 5 - 10%. 
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3.2 Effective range corrections 



In next order of the momentum expansion of the effective field theory, there is no operator 
which induces S — D mixing of the deuteron state. It will first appear at one order 
higher The dimension-eight couphngs C2 and C| give the additional diagrams shown 




Figure 8: Corrections to the fusion amplitude coming in at next-to-leading order, 
in Fig. 8 in first order perturbation theory. Each such operator V2 has a momentum 
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matrix element (k|V2| q) = C2(k^ + q^)/2. The contribution from Fig. 8a is seen to be 



This can be expressed in terms of the divergent integral 

r H^k 1 —M 

which is the same as occured in the lowest-order determination of the wavefunction renor- 
malization constant in (Q). It is finite after PDS regularization which gives for the other 
occuring integral 

^2(7) = -mJ (03^^ = -7^/0(7) (36) 

since / d'^k/{2TT)'^ = in dimensional regularization. Together with the function Aq^p) in 
(^) and the related function 

we thus have for the matrix element ( p^ ) 

Ta = \ci[h{i)Ao{p)+h{l)A2{p)] 
In the same way as we could express the integral 12(7) in terms of loij), we also find 

A2{p) = MiPp{0)-j^Ao{p) (38) 

Here ipp{0) = CrjC^'^^ where ctq is the Coulomb S-wave phaseshift. When we eventually 
use this result to calculate the fusion rate from (||), we will take the absolute value and 
this phase factor will not contribute. We therefore write 

Ta = Ic^Ioh) [MCrjip) - 27'^o(p) 

where the same phase factor also should be dropped in the last term. This result is now 
to be taken in the fusion limit p ^ as in the previous section. 

The contribution from the diagram Fig. 8b involves the Coulomb Green's function and 
its derivative in the triple integral 



resulting from just one proton bubble. This can be expressed in terms of the function 
Bq{p) in (24) and the related function 

r d^k d^n k^ 
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as 

ni = \ctci[h{i)B^{p) + h{i)B2{p)]M'^) 

With the simphfication 

B2{p) = MJo{p)-^^Bo{p) (40) 
we find the total contribution 

= ^^' l-^^Joi) 1^"^°^^^ " 27^^oW] C, (41) 

from all the Coulomb-dressed proton bubble diagrams in Fig. 8b. Here we have again 
replaced V'p(O) by Cr^{p). 

The remaining diagrams involve the proton derivative coupling Diagram Fig. 8c 
gives 

This can again be expressed in terms of the functions Bq{p) and Cr;(p) and their derivatives. 
In particular, we define 

and introduce 

Mp) = J (^kVp(k) (44) 

which is the double derivative of the Coulomb wavefunction at the origin. Both of them are 

highly divergent, but can be calculated in the PDS regularization scheme and expressed 
in terms of already introduced functions. This is shown in the appendix. We thus find for 
this diagram 

Tc = \cl [Mp)b'2{p) + Mp)Bo{p)] 

where ip2ip) = —aMiJ,'ipQ{p) in the limit p — > as shown in the appendix. Then we also 
have B2{p) = MIq{^) — aMiJ,BQ{p) and we therefore get 

Tc = [M/o(7) - 2aM/x5o(p)] (45) 

which now involves only finite and known quantities. 
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In Fig.Sd we sum over all the Coulomb-dressed proton bubbles. The result is given 
by the multiple integral 



M 



X Gc(i^;k',k)j-^2^^p(0) 

Again we can reorder the integrand so that the result is expressed in terms of simpler 
functions, 

where now 

involves the derivative of the Coulomb propagator. It is also evaluated in the appendix. 
In the limit p — > we find J2{p) = —aMnJ^ij)) which together with the related result for 
B'^ip) gives 

T, = \CI^^^^^^ WU^) - '^o.M^,B,{p)] (47) 

It has the same structure as Tc in (^) and they can therefore be combined into a simpler 
result. 

The deuteron side of the diagrams in Fig.Se is seen to be just Bq{p). Summing up the 
bubbles on the proton side, we find 

Te = ll^iipj ['^2(P)V'0(P) + Mp)Mp)] BoiP) 
In the limit p — > this simplifies again with the result 

Te = -CPCr,j^^^l^aMf,Bo{p) (48) 
Similarly we find that the diagrams in Fig.Sf gives 

= (^ i_cpj,{p) ) ['^2(p)^o(p) + Up)J2(j})] Mp)Bo{p) 

which becomes 



cIMp) 



2 



= -CICAyt^j^A aMfiBoip) (49) 
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in the low-energy limit p ^ 0. 

Adding now up the contributions from all the diagrams in Fig. 8, we obtain the sum 



1 



. MIoij)Cr, 



(1 - Co%(p))2 



o'^o(p) 



where Tq{p) is the lowest order matrix element (p5| 
to this order is therefore 



The full transition matrix element 



Tfi = Z2To{p) + ZoT2{p) 

where Z2 is the next-to-leading order renormalization constant ( [lO| 
combining terms, we obtain 



^0 ^Tfi 



+ 



Aoip) + Bo{p)Cr, 



aMfiC^ 



1-CIUp) [i-cIUpW 
Ao{p) + Bo{p) 



1 - c^Mp) 



Reordering and 

(50) 
(51) 
(52) 



In the bubble integral Jo{p) we can take p — > since it is finite. The function Bq{p) 
is also finite in this limit while Aq{p) becomes proportional to the Coulomb factor Crjip) 
which diverges. As shown previously in the application of the same effective theory to low- 
energy, elastic proton-proton scattering, the first square bracket is now just the physical 
proton-proton scattering length Op calculated in next-to-leading order with the result |lj 



M 
4tt 



aMixCl 



'^0 



(53) 



The last term is the effective-range correction which is important in order to have a 
physically meaningful result for the scattering length. We see that when this is zero, we 
have the previous result (15) used in leading order. 

The transition matrix element in next-to- leading order is now given by (^2|). Isolating 
a common factor, the reduced matrix element (^) follows as 



A2(0) = Ao(0) 



Cf(/x) + C|(A.) 



(54) 



where Ao(0) is the leading-order result ( pSD but now expressed in terms of the next-to- 
leading order scattering length (|5^). With the already established value for the coupling 
constant C2 we see that it is now multiplied by the factor 1-1-7^^/2. This can be interpreted 
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as the first term in the expansion of the deuteron normalization factor discussed in the 
previous section. While this term in the result is independent of the renormalization scale 
/i, we see that the last term is generally not. However, when [i is much larger than the 
other mass scales given by the scattering lengths, this dependence goes away and we are 
left with the definite result 



A2(0) 



M>7 



Ao(0) 



1 

1 + ipPd 



(55) 



It has a structure which is very similar to the reduced matrix element in the standard 
nuclear physics effective-range approximation [0] [0] 



Ao(0) + ^ap72(/9p + prf) 



(56) 



With the known values for the different nucleon parameters, the result in this old ap- 
proximation is therefore A£;/{(0) = 2.66 or A|;^(0) = 7.08. On the other hand, our 
next-to-leading order result (|55|) gives A2(0)^>^ = 2.54 which is just a 1.4% addition to 
the leading order result we previously obtained . This is surprisingly small, but results 
from an almost total cancellation between the two effective-range corrections in (|55|). The 
net result for the squared matrix element is A2(0)^>^ = 6.45 which is seen to be 8% below 
the effective-range value. 



3.3 Contribution from counterterm 



A complete calculation of the fusion rate in next-to-leading order must include all operators 
contributing to this order in the momentum expansion of effective theory. Until now we 
have only included the effects of the dimension-eight operators coupling four nucleons with 
a derivative interaction. Since our result above in general depends on the renormalization 
scale, it signals that the calculation is incomplete. There should be additional interaction 
terms that in principle should absorb all dependence on the renormalization scale. This 
is in fact the shown by Butler and Chen[pO| and discussed in the introductory 

section. It has the structure as given in ( p^ ) and corresponds to the weak current coupling 
directly to the four-nucleon vertex. In a more fundamental theory it could be due to 
weak interactions via virtual pions, coupling to excited nucleons or more general two- 
body operators in nuclear physics language. Obviously, this counterterm will also modify 
the numerical result for the fusion rate in addition to softening the //-dependence. 

In our case it gives a contribution depicted by the Feynman diagram in Fig. 9a. It is 
similar to the previously calculated contribution from Fig. 8a in (34) and becomes 

(2^ J (4^i?T^^-(^) = -^^^^M^^ (57) 

The strong interactions in the initial state, now to lowest order in the derivative expansion, 
gives the series of diagrams shown in Fig. 9b. They form again an infinite geometric series 
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Figure 9: Contributions to tlie fusion amplitude from the counterterm. 



whose sum 



Ti ct 
h 



Li A C, 



cIMp) 



(58) 



cIMp) 

is given by the proton-proton physical scattering length Op from ( [T^ ) and a(/i) from ( pO| ) 
in the fusion limit p ^ 0. With the regularized value for the integral Iq, we thus find the 
total contribution from the counterterms to be 

1 



T ct 



T r ^( 



7) 



(59) 



The corresponding reduced matrix element then follows from (^) after multiplication by 
the wavefunction renormalization constant Zq. 

We now include this new contribution as a correction to the matrix element ( |5^ ) 
coming from the ordinary axial current interactions. For the combined result we then 
have 



AS* 



(0) 



Ao(0) 
Lia(/u) 



1 

1 + 2^Pd 



(/^ -t) 



1 



a(/x) 



^(C7f(/.)+C2^(/.) 



(60) 



The coupling constant Lia of the counterterm must have a dependence on the renormal- 
ization scale iJ, so that the total //-dependence in the last term is negligible. When 3> 7 
we see that this requirement leads to 



(61) 



where d-iA is an unknown, dimensionless constant. It is set by physics on scales shorter 
than l/rriT^ and its natural value should be around one as pointed out by Butler and 



Chen|2C]. In order to get a rough idea of the sensitivity of the result on this parameter, 
we take ^ = m.,^ which is the scale at which one should match the effective theory to the 
more fundamental theory involving pions. Varying then Iia in the interval [—1, 1], we find 
that the fusion rate measured by A2*(0)^ varies linearly from 6.22 to 6.84. These values 



are seen to be systematically below the effective-range result following from (|56|), but are 
within the 5% - 10% uncertainty range expected at this order. 
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4 Discussion and conclusion 



Effective field theory is a very powerful approach to low-energy phsyics. It can hardly 
be said to be wrong when used correctly since it is just based upon the basic symmetries 
of the problem and standard quantum field theory. In that way it is a very conservative 
approach since it does not admit assumptions about the physics on scales shorter than 
it is meant to handle. Instead of such specific and model-dependent assumptions, one 
has higher-dimensional contact interactions and counterterms with coupling constants 
which represent the unknown physics. The most common criticism against effective field 
theory is therefore that it is not accurate enough since the results may depend on one 
or more such coupling constants which are not a priori known. One can make estimates 
of these unknown coupling constants based upon some kind of naturalness supported by 
dimensional analysis and the renormalization group. 

But these counterterms do not really represent a weakness of effective field theory. 
Since they are interactions appearing in a Lagrangian, they will appear with the same 
strength in many different processes. If one or more of these allow for the determination of 
the corresponding coupling constants, one can then make much more accurate predictions 
for the other reactions. One recent example is radiative neutron-proton capture n + p ^ 
d + J. When the process takes place at very low energies or at rest, it is dominated by 
a magnetic dipole transition which at next-to-leading order also involves a four-nucleon 
counterterm very similar to the one we have considered here for proton-proton fusion. 
From the measured rate at these low energies, the counterterm can then be determined 
numerically [p!^] . The same neutron-proton fusion process is also a key reaction in big-bang 
nucleosynthesis where it takes place at energies upto around 1 MeV. Chen and Savage have 
now calculated the corresponding cross section with an uncertainty of 4% based on the 
measured counterterm |]2^ . A similar accuracy can be expected also for proton-proton 
fusion if the counterterm can be determined in some other process. 

It has already been pointed out that our results for the proton-proton fusion have a 
very similar structure to what one finds in the effective-range approximation in nuclear 
physics. This has also been seen in other processes investigated within the same effective 
theory and at higher orders in the perturbative expansion [|T2|| . It is understood when 
one realizes that these processes are dominated by the properties of the deuteron wave- 
function at large distance scales which is contained in the effective-range approximation. 
In the KSW field theory, these properties are coded into the coupling constants Cq and C2 • 
While Cq is responsible for binding the deuteron and must be treated non-perturbatively, 
the effects of to be treated perturbatively and gives the detailed behaviour of the 

wavefunction at large distances. In the above C2 was determined by matching to the effec- 
tive range parameter pd- In order to get better agreement with low-energy proton- neutron 
scattering data which are related directly to the deuteron bound state wavefunction via 
analytical continuation, it has recently been pointed out by Phillips, Rupak and Savage 
that one should instead match to the wavefunction normalization parameter Z^ pl]] . 
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This gives the result 



where = 1.69. They have shown that this markedly improves the convergence of the 
perturbative calculation of many processes involving deuterons at low energies. Rupak has 
recently applied this improved method to neutron-proton fusion n + p ^ d + 'j at energies 
relevant to big-bang nucleosynthesis as discussed abovepg]]. Including one higher order in 
the perturbative expansion of the elctric transition amplitude, he has then obtained an 
accuracy of 1% for the calculated cross section. 



In our case we can now use this new value for C2 in the result ( |5^ for the reduced 
matrix element. Including also the counterterm as in (|60|), we then obtain our final result. 
Again the counterterm coupling constant will have the form (|6l| ) for large values of the 
renormalization mass. Choosing = itIt^, we now find that A|(0) varies between 7.04 and 
7.70 when the parameter iiA takes values in the interval [—1,1]. With the size of the 
unknown counterterm in this range, we thus have the central value A|(0) = 7.37 with a 
conservative estimate for the uncertainty of 6% - 8%. We thus find a somewhat higher 
value for the fusion rate in this improved perturbative calculation compared with results 
from effective range theory ( [5^ ) and the inclusion of axial two-body effects |6|. But these 
other results are now within the accuracy range of the effective field theory method. Should 
there in the future turn out to be a real discrepancy between these different theoretical 
descriptions, it can only be due to meson processes at shorter scales which have been 
overlooked or not correctly handled in the more conventional nuclear physics approach. 

The only way to improve the accuracy and thus obtain a more predictive result, is 
to establish the value of the counterterm coupling constant Lia- In principle it could be 
measured in many other reactions, but the most promising is elastic and inelastic scattering 



of neutrinos on deuterons as shown by Butler and Chen|2C]. High-precision experimental 
results for these reactions would then result in a known value of the relevant counterterm. 
As shown by Rupak for radiative neutron capture, one should then be in the position to 
obtain the proton-proton fusion rate with a much improved accuracy. This will place our 
understanding of this fundamental process on a more solid basis. Needless to say, it will 
also strengthen our knowledge of the neutrino production rate in the Sun. 
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6 Appendix 



We will here regularize and evaluate the divergent integrals involving Coulomb wavefunc- 
tions which are needed for the effective-range corrections to the fusion rate. Some of 
them have previously been encountered in connection with higher order corrections to 
low-energy proton-proton elastic scattering[|l7|]. They were then calculated by a method 
based on regularization of the Fourier-transformed Coulomb wavefunctions. We will here 
use a different and simpler method. 

The simplest integral is J2{p) in (|4^) which we rewrite as 

J2(p)=/Jo(p)+|^/^(k2-p2)(k|Gc(i?)|q) 

It represents a Coulomb-dressed bubble propagator with a derivative interaction at one ver- 
tex. Here we have introduced the free eigen-momentum states (k | and | q) . The Coulomb 
propagator Gc{E) satifies the Lippmann-Schwinger equation Gc = Gq-\- GqVcGc where 
Go{E) is the free propagator ( p!^ ) and Vc is the Coulomb potential. In momentum space 
it has the matrix element (klVclk') = 47ra/(k — k')^. The first term will now give zero 
with the use of dimensional regularization, 

We then insert two complete sets of momentum eigenstates between the three operators 
in the matrix elements in the second term. The denominator in the free propagator Go 
then cancels against the factor k^ — in the integral. We are thus left with 

In the integral over k we now shift the integration variable k ^ k + k' and use the PDS 
regularization result 

d^k ina 



(27r)3 k2 



an (64) 



The remaining two integrals over k' and q then simply gives Jo{p)- We thus have the 
result 

J2{p) = \p^-aMfj]Joip) (65) 
Except for a higher order term in the fine-structure constant a, this agrees with what we 



obtained with the much more cumbersome wavefunction regularization method[17|. 



The next integral ip2{p) in (0) corresponds to the double drivative of the Coulomb 
wavefunction at the origin. We can write it as 

d^k 



^(k2-p2)(k|V'p 
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where | V'p) is a Coulomb state with momentum p. It can formally be expressed in terms 
of the free state | p) as 



\i^p) = [l + GcVc]\p) 



One then has 

^(k2-p2)(k|V^p) = J^(],^-p^)[{],\GoVc\p) + {k\GoVcGcVc\p)] 

= /(03(k^-P^)(k|GoyclV'p) 

using GcVc\p) = \ ipp) — \ p) in the last term. Inserting now again two complete sets of 
free momentum states as above, it follows that 

mp) = P mp) -mJj^j (2^)3 (k-q)2 ^q'^p) 

After a shift of integration variable, we have the result 

Mp) = \P^ -aMfj]i>o{p) (66) 
when making use the the PDS regularized integral (|6^. 

The last integral we need is B'2{p) in (|43|). Rewriting it as above, it takes the form 

where the first term is the finite integral (24). In the second term we can use the Lippmann- 
Schwinger equation for the Coulomb propagator. Again we find that the denominator of 
the free propagator cancels against — p^ in the numerator. The first term in the integral 
gives then just the integral /o(p) in (pSj). Going through the same steps as above with 
insertion of complete sets of states, the second term is then reduced to the finite integral 
Bq{jp). In this way we obtain 

B'^ (p) = MIo (p) + [p2 _ aM^] Bo (p) (67) 

which again is a surprising simple result. 

We notice that these three divergent Coulomb integrals contain the common factor 
p^ — aM^jL in the results. This can be understood as coming from the divergence of 
the double derivative of the Coulomb wavefunction ■0p(r) at the origin. It satisfies the 
Schodinger wave equation 



where the energy E = p /M. When we now take the limit r ^ 0, it follows that 

- V2V;p(r),^o = [P^ - aMfi] V^p(O) 
since the regularized integral (|64[) is just the Coulomb potential at the origin. 
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